clc;

a = 0.35;

L = 1.0;

c = 4.0*pi*a*a*a/3/L/L/L;

poro = 1.0 - c;

cmax = pi/6;

kai = (c/cmax)^(1.0/3);

% FCC
% q = [0.1, 0.1620994, 0.2627620, 0.3518875, 0.4503759, 0.5354862, 0.6240194, 0.7048893, 0.7778734, 0.8380856, 0.9093106,...
%      1.004412, 1.099079, 1.176754, 1.234515, 1.261369, 1.271545, 1.295785, 1.297464, 1.259288, 1.233984, 1.238161, 1.310988, 1.251041, ...
%      1.089836, 1.059025, 1.210762, 1.221453, 0.9163566, 0.5120200, 0.4523067];

q = [0.1, 0.1418649, 0.2012564, 0.23331523, 0.2564809, 0.2584787, 0.2873609, 0.3340163, 0.3536763, 0.3253622, 0.2689757,...
     0.2037769, 0.1809341, 0.1877347, 0.1534685, 0.09034708, 0.02857896, -0.05512626, -0.1278724, 0.1013350, 0.5492491, 0.4615388, -0.05736023, -0.2865924, ...
     -0.4709215, -0.6870076, 0.1455304, 1.251891, 0.9742811, -0.8566269];
 
 q = 10.0*q;
 
 len = length(q);
 
 Ds = 0.0;
 
 for i = 1 : len;
     
     Ds = Ds + q(i)*kai^(i-1);
     
 end
 
 Ds
 
 d = 2.0*a;
 
 K = L*L*L/6/pi/a/Ds